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ABSTRACT 



In 2D-simulations of self-gravitating gaseous discs, the potential is often computed in the framework of "softened gravity" initially 
designed for A'-body codes. In this special context, the role of the softening length A is twofold: i) to avoid numerical singularities in 
the integral representation of the potential (i.e., arising when the separation |r - r'| ^ 0), and ii) to acount for stratification of matter 
in the direction perpendicular to the disc mid-plane. So far, most studies have considered /I as a free parameter and various values or 
formulae have been proposed without much mathematical justification. In this paper, we demonstrate by means of a rigorous calculus 
that it is possible to define A such that the gravitational potential of a flat disc coincides at order zero with that of a geometically thin 
disc of the same surface density. Our prescription for A, valid in the local, axisymmetric limit, has the required properties i) and ii). It 
is mainly an analytical function of the radius and disc thickness, and is sensitive to the vertical stratification. For mass density profiles 
considered (namely, profiles expandable over even powers of the altitude), we find that /I : i) is independant of the numerical mesh, 
ii) is always a fraction of the local thickness H, iii) goes through a minimum at the singularity (i.e., at null separation), and iv) is such 
that 0.13 < A/H < 0.29 typically (depending on the separation and on density profile). These results should help us to improve the 
quality of 2D- and 3D-simulations of gaseous discs in several respects (physical realism, accuracy, and computing time). 

Key words. Accretion, accretion discs — Gravitation — Methods: analytical — Methods: numerical 



1. Introduction 

The computation of gravitational potentials and forces is a criti- 
cal step in many astrophysical problems where gravity plays an 
important role. Because Newton's law of gravitation diverges as 
the relative separation |r, - rjl between particles vanishes, one 
classically adds a positive constant A referred as the "soften- 
ing length". This technique of singularity avoidance, initially 
developped to prevent binary collisions and particle evapora- 
tion in A^-body simulations, has led to the c oncept of "soft- 
ened gravi ty" (Hocknev & Eastwood 1988; Ant hony & Carlber^ 
ll988inSommer-Larsen et al.lll998l: INelsonll2006 ). Soon, people 
realized that the use of a softening length introduces a notice- 
able bias in models, especially for the stability properties of 
stellar systems, and there have been several attempts to search 
for the most appropriate /t- value (e.g., Romeo 1994, 1997; 
ISommer-Larsen et al Jl 1 99i lDehnenll2()0 ll) . 

Softened gravity has also been employed in the compu- 
tation of the gravitational potential of gaseous discs (e .g., 
Papaloizou & Lin 1989; Adams et al. 1989; Saio & Yoshii 1990; 
Shu et al. 1990; Morishima & Saio 1994; Sterzik et al. 1995£ 



Laughlin et al III 9971 Il998t lTremaiiiJl200lt ICaunt & TaggS 



20011: iBaruteau & Massetl 120081; iLi et af] l2009t) with the addi 



tional justification that the softening length takes into account 
the vertical stratification of matter. Various prescriptions, gener- 
ally in the form of a function of the cylindrical radius and/or of 
the disc parameters, have been adopted (see references above) 
without convergence towards a "universal formula". As for sys- 
tems of particles, it is recognized that the softening length can 
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dramatically affect the stability of gaseous discs (again, see ref- 
erences hereabove). It is therefore fundamental to ask i) whether 
or not softened gravity can definitively help to determine or 
mimic the Newtonian potential of a gaseous disc with a cer- 
tain level of accuracy, and ii) eventually which formula — if 
one exists — is the most appropriate. This is the aim of this pa- 
per. By equating the softened potential of a flat disc to that of a 
geometrically thin disc, we find that we can define a softening 
length A which has the required properties. We therefore report 
the first reliable formula for the softening length on the basis of 
approximate but rigorous calculus. In the axisymmetric limit, A 
is found to be a sharp function of the relative separation \r - r'\; 
it is a fraction (of the order of 4: at the singularity) of the disc lo- 
cal thickness (see also Hure & Pierens 2006; Baruteau & MasseJ 
200 81), and it does not depend on the numerical resolution (see 
Li et al. I l2009h . The formula, easy to implement into hydrody- 
namical codes of self-gravitating 2D- and 3D-discs, will enable 
to increase the degree of reaslim of simulations, both by preserv- 
ing the Newtonian character of the potential and force field, and 
by accounting for the vertical structure. 

The paper is organized as follows. In Sect. |2] we recall the 
expression for the midplane gravitational potential of a geomet- 
rically thin disc as well as that of a flat (i.e., zero thickness) disc. 
Because of the hypothesis of axial symmetry, the integral ker- 
nel contains a complete elliptic integral of the first kind that can 
be appropriately expanded at the singularity |r - r'| - and 
its neighborhood. In Sect. |3] we estimate the effect of vertical 
stratification by performing the analytical integration along the 
z-direction. Various mass density profiles corresponding to fi- 
nite size discs are considered, namely the homogeneous case and 
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mixtures of even powers of the altitude. In Sect.|4l we show that 
this calculus naturally leads to a local prescription for the soften- 
ing length A. We present a table showing the great diversity (and 
incoherence) of prescriptions used so far and derive the general 
relation for A. We discuss the sensitivity of the softening length 
to the separation and disc aspect ratio. In Sect.|5] we summarize 
the results and suggest a few possible extensions and generaliza- 
tions of this work. 



2. Mid-plane gravitational potential in both flat discs 
and geometrically thin discs : local treatment. 

2.1. Notation and background 

We consider two axially symmetric discs (see Fig. [TJ : i) a flat 
(i.e., zero thickness) disc with inner edge Oin, outer edge Oout, and 
total surface density Stot; and ii) a geometrically thin disc with 
the same edges and same total surface density, but local thick- 
ness H - 2h > 0. For the geometrically thin disc, the condition 



(If 



«: 1 



(1) 



is assumed at any radius a e [flin, flout] (lPringIell98lh . Moreover, 
we have 

p(z)dz = Stot, (2) 



r 



where p is the mass density at the altitude z from the mid-plane, 
Z- is the altitude of the bottom of the disc, and z+ = z-+ 2h is for 
the top. In general, p, Stot, Z-, and h are expected to depend on 
the radius a. In the following, we shall also consider symmetry 
with respect to the mid-plane, which is the case in most disc 
models (not true for warped discs for instance). This means that 
Z- + z+ = 0, and z+ - h. For z-i- — > for all a, the two discs 
therefore become equivalent. 



a geometrically thin disc 



a,R 




a flat disc 



Fig. 1. A flat disc and a geometrically thin disc, axially symmet- 
ric and symmetric with respect to the mid-plane. The disc edges 
are flin and flout, h is the semi-thickness, Stot is the total surface 
density, and a and R are cylindrical radii. 



For the thin disc, the mid-plane gravit ational potenti al at ra- 
dius R is given by the double integral (e.g. lDuran dl[T95l 



*P*in(^-) ^ _2(5 pkK(k)dz, (3) 

Jfli,, V R J-h 

2 



where 



k = 



V(fl + R)^ + z^ 



(4) 



is the modulus and K is the complete elliptic integral of the first 
kind (Legendre form). For the flat disc, we have p{z) = S,ot5(z) 
everywhere. Consequently, the potential at radius R is 



xpflat^^) 



Stot'nK(m)Jfl, 



where 



2 

a + R 



(5) 



(6) 



is the associated modulus (note that k - m for z - 0). Both Eqs. 
|3] and |5] involve a logarithmic singularity when the modulus of 
the elliptic integral approaches unity. In spite of this divergence, 
the potential is generally finite at any radius. It is interesting to 
see that Eq. |3]can also be written as 



Vj/thin^-^-j 



-2G 



where 



^^totA" 



J-h 



pkK(k)dz. 



(7) 



(8) 



Clearly, is a function of a, R, and h and contains the effects 
of vertical stratification. So, if x is known preliminarily at all 
radii and for a given disc structure (edges, thickness, mass den- 
sity profile), ^F*'" is found from a single integral over the ra- 
dius (as for *!'*''"). It is therefore advantageous to have a one- 
dimensional formula that describes a bi-dimensional distribu- 
tion of matter (in the axially symmetric case), especially in terms 
of computing time. By comparing Eqs.|5]and|7l it is tempting to 
setx - mj^im^), where m^ is a certain modulus to be determined 
(see Sect.©. 

2.2. Expansion around tlie singuiarity 

The presence of the function K does not allow us to derive x 
analytically for any mass density distribution. The expansion of 
K(^) appropriate for treatin g the logarithmic singularity is (e.g. 
iGradshteyn & Ryzhiklll965h 



K(^) = YjP„Pn{k')k'^ 



(9) 



«=0 



where k' - Vl - k^ is the complementary modulus, and 



Po = ln^, 



Pn^l(k') = P„(k') - 



I 



(2«+l)(n+l) ' 



for n > 0, 



(10) 



By construction, this expansion is efficient when k' — » 0, which 
corresponds to A; — » 1 . It can be shown that the condition k'^ «: I 
implies that 



h 



and 



a -¥ R 

\R-a 



<S 1, 



fl 

«: -. 
h h 



(11) 
(12) 
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The first inequality is automatically fulfilled within the geometri- 
cally thin disc approximation (see Eq.[TJ- The second inequality 
means that the present calculus is valid only at the singularity 
R — a and /or a few disc thicknesses in radius. This is precisely 
the radial d omain where the numer ical determination of ^F*'" is 
tricky (e.g.. IStemwedel et alJll990l) . because of the divergence 
of K(A;). We note that there is no special constraint on the local 
shape h{a) of the disc provided that this remains geometrically 
thin. 

To the lowest order, K(A;) ~ In p- and A: » 1, and Eq. [8] be- 
comes : 



with an error of 0{k' \nk'). 



r 4 

pln—dz, 
Jo k' 



(13) 



3. Effect of vertical stratification : an estimate of tlie 
;(^-f unction 

3. 1 . Vertically homogeneous discs 

We can estimate x from EqlT3]in a few particular cases, for in- 
stance when the disc is vertically homogeneous (x is denoted ;\fo 
in this case). If p does not depend on z, then Etot = ^ph = Eo and 
we find that 

Xo = ln4 - In^^ - In/o, 



where 



k'- 



(14) 



(15) 



In/o - ^arctan/7_ - ^arctan;7+. 

This is a complicated function of R/a, where the local aspect 
ratio 2/j/a is the unique parameter If we define the variable 



R-a 



which measures the radial separation from the singularity R-a 
in units of the disc semi-thickness h, and 



e = 



2a' 



(17) 



which is the quarter of the aspect ratio, we find that 
XQ = ln4-ln 



(1 -Hx2)(l -Hejc)2 



1 1 + ex e 

- X arctan — i arctan . 

X e \ + ex 



(18) 



The function ;if(x) is plotted in Fig.|2]for hfa - 0.1 (that is e 
0.05). To order zero around x = 0, Eq. [I4lreads 



4 n I 

1 +\n- + ex \x\ - In Vl -H x^ -H jc arctan x. (19) 

£ 2 



We see that xo reaches a maximum value of 1 + In ^ » 1 at 
R = a, and is slightly asymmetric with respect to x = 0, since 
we havexoix - +1) ~ Xo(x - -1) + 2e. 




-3 -2-10 1 

X 

Fig. 2. Variation ofxo with x according to Eq.fT4l 



3.2. Effect of vertical stratification 

We can probe the sensitivity ofx to vertical stratification by con- 
sidering the following profiles : 



2q 



for |z| < h, 
otherwise. 



(20) 



where po is the density at the disc mid-plane (possibly a func- 
tion of a) and ^ > 1 is an integer. These correspond to finite size 
discs. When q - I, the profile is close to a Gaussian distribution, 
a case found for instan ce in vertically isothermal disc at hydro- 
static equilibrium (e.g.. Chi ang & Goldreichiri997l; iHirose et al.l 
2006; Edgar & Guillen 2008). As q increases, the vertical profile 
becomes flatter and flatter For ^ — > co, we recover the homoge- 
neous case considered above. 

We can calculate again from Eq.[T3] but using Eq.|20l From 



(16) Eq.m we get Etot = ^'^o, and then 



2q+\ 



(21) 



where we have defined 



Xq 



= (2^+1) j| g)''^K(^)«'|. (22) 



We can perform the integration using the expansion of K(A;) 
around k' - 0, as above. To the lowest-order, we have 



Xq = \nA-\\\k'^-\\\fq, 



where 



In/, = {-If 



arctan?/- arctan/y^ 



2o+I 



2q+l 



(23) 



(24) 



2(k-q) 
k'l- 



2(k-q) 



k=0 



2k +1 



This expression can be rewritten in a different form since 
each sum represent the first q + I terms of the expansion of the 
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arctan function. In particular, a form appropriate to numerical ^ ^ 
computations around the singularity is 



2,1 ^ 11 

arctan — 



(~1)* 2(l<-q) 



(25) 



+ 77: 



k=0 



k„2k 



2{k + q) + 3 



We see that;^f^ is a function of R/a, and h/a is the parameter. 
Using again the variable x and the e-parameter, we find 



Xq =ln4-ln 

-{-If 



(1 +x^)(l + ex)^ 



+ (I + ex)2 

HI' 



(26) 



X arctan x 



\l+ex) ^^2(k + q) + 3\l+ex) ' 

The function is displayed versus x in Fig.[3]for h/a = 0.1 and 
different values of q. Around x = 0, Eq.l26lbecomes 



1 



Xq 



-I- In - + ex - In Vl + x^ 
2q+\ e 



and we note that 
2q 



Xo-Xq 



2q + 



-(-l)V«(^|x|-xarctanx), (27) 
- + [(-l)V'? - 1] - X arctan . (28) 



The maximum of Xq still occurs at ^ = a (see Fig. [3]) 

2q 



2q+l' 



(29) 



Again, we note that the difference;t',(/? = a+h)-Xq{R - a-h) ^ 
2e is insensitive to q at the actual order of precision. 

Finally, for the vertical profile defined by Eq. |20l x follows 
from Eqs.[T8]|2Tl and|26l We see that;^' reaches a maximum at 
R - a, and 

;rmax * 1 + In - + (30) 
e 2q + 1 

for q > 1. For large values of the ^-parameter, ^ Xo- Figure 
|4] displays x versus x for h/a - 0.1 and different values of q. 
It is remarkable that x is very weakly sensitive to the vertical 
profile and that x remains very close to xo (the homogeneous 
case), especially for |x| > 1 . 



3.3. Full generalization 

The generalization of the above result is relatively straightfor- 
ward for profiles of the form 



CO 

p(z)=Po2'^«(^) ' 

n-0 



(31) 



where c„ are real coefficients. In this case, we find that 



X = 



y±!^, withEtot = 2;oy — ^, (32) 





1 




1 ' 1 ' 1 ' 1 
q=0 (homogeneous case) 


- 






/ \ 1 (parabolic) 
/ /\ 2 




~h/a 




/v^riS<v 3 

/|>2^^^v^ - 






'Afy^ CO 
















■ ~h/a 






R=a-h 


II ^ 




r , 1 


1 , ,1,1 



Fig.3. Variation of;^'^, with x according to Eq.|23]for h/a - 0.1 
and different values of q, namely ^ = (the homogeneous case), 
q- \ (the parabolic case), ^ = 2, 3, 5 and ^ — > 00. 



-5.715 
-5.382 



-5.753 



1 (parabolic) 
Cosine 
2 




Fig. 4. X versus x for the mass density profile defined by Eq. 
The cosine profile discussed in Sect |4.2| is also shown. 



where Xn is found from Eq. |26] Again, is a function of x, and 
e is the parameter. It is especially convenient to write x as 



X =^0 + 5^0, 



(33) 



where 5x(^ is the deviation relative to the homogeneous case. We 
note that dxo is generally small compared to x^ ^s long as the 
bulk of the mass is localized in the disc mid-plane. We conclude 
that;^' is determined mainly by the homogeneous contribution xa, 
provided that the matter is gathered around the mid-plane. For 
X — > 0, we get from Eqs.[T9]andl27] 



^ Cn { 2n 



^° ~ ^ 2n -H 1 I 2n + 1 

«=() 



(34) 



[(-l)"x-" - l]^^|x| -X arctan xj 

We see i]\2ix reaches a maximum value at x = 0. This value is 

2nc„ 



2o 

2, 



. {2n + 



„ot — (2«+l)2- 



(35) 
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4. Application to softened gravity : a prescription 
for the softening length 

4.1. The avoidance of point mass singularities 

As mentioned, the complete elliptic integral K exhibits a log- 
arithmic divergence as its modulus — > 1. This is the direct 
consequence of the point mass singularity (i.e., |r - r'\ — > 0). 
The determination of accurate potentials from Eqs. |3]and|5]is 
therefore not straightforwa rd and requires a careful treatment 
of improper integrals (e.g., IStemwedel et al.l[T990l: lHurell2005t 
iHure & Pierensll2005b . This technical difficulty is usually cir- 
cumvented by changing the relative separation according to 



Reference 



softening length A 



(36) 



where /i is a constant known as the "softening length". 
The main drawback is that softened gravity modifies Newton's 
law for gravitation both on short and long ranges. It lowers the 
magnitude of forces, enhances stability, and introduces a bias in 
models that is not easy to measure an d interpret (for stellar and 
gas discs see, e.g., Papaloiz ou & Lin|[l 989: Saio & Yoshii 199(J 
lRomeolll994l: [SCTnmer-Lars en et alJll998ciDehnenll2001,) . 

In the case of gaseous discs of interest here, the modifica- 
tion of the relative distances according to Eq.[36]changes the ex- 
pressions for ^i*"^' and ^F^m. ft is however easy to show that the 
associated "softened" potentials denoted ^F^"' and ^F*'", respec- 
tively, can still be written in terms of Eqs. [5] and |5j respectively, 
provided that the modulus k is replaced by 



and m is replaced by 



2 



2 



V(fl + R)^ + 



(37) 



(38) 



In disc models and simulations, one never (or rarely) computes 
(jjthi n^ or its fully asymmetric/tri-dimensional version (see how- 
ever lLi e t al. 2009). Instead, one computes ^I*"^* in the framework 
of softened gravity, that is ^F^^*. The softening length A appear- 
ing in Eq.[38]must therefore be prescribed. We note that solving 
Eq.[38]for/lleadsto 



4aR (a - R)^ 



(39) 



where m'^ - yj\ — ml is the complementary modulus. 

4.2. A prescription for the softening iength A 

Many prescriptions for A have been proposed. In general, this 
is not a constant but a certain function of the radius and/or disc 
parameters. Table [T] gathers a few formulae for A used by dif- 
ferent authors over twenty years. Although not exhaustive, this 
list clearly shows that there is no trend in magnitude and varia- 
tion in space (and possible dependency on the disc parameters). 
The results obtained in Sect. [3] can help substantially to define 
the appropriate prescription for A. We can see from Eqs. |5] and 
|7]that the gravitational potential in the mid-plane of a thin disc 
is equivalent to the softened potential caused by aflat disc pro- 
vided that A is the root of the equation 



pkK(k)dz - 0, 

h 



(40) 



Papaloizou & Lin (1989) 
Adams et al. (1989) 
Saio & Yoshii (1990) 
Shuetal. (1990) 
Morishima & Saio (1994 ) 
Sterzik et al. (1995)° 
_Laughlin et al. (1997) 
Laughlin et al. (1998) 
Tremaine (20(31) 
Caunt & Tagger (2001) '' 
Baruteau & Masset (2008) 
Li_e t al. (2009)"^ 

this work'' 



0.056 X (Oout - am) X f(a,R, a,„,aoM) 
OAxR 

(0.01,0.02) XfloLit 
OAxR 
0.02 X Ooi,, 
0.06 X A^ 

0.1 xai„,0.1 X ao„„ and 0.001 xR 
0.0001 + O.OIR X f(a, a,„, 
/3xR, withySa; lO"'* - 0.2 
H = 2h 

0.3 - 0.5/i (depending on scale height) 
«O17to0.33xAa 

X h/e alR = a, homogeneous case; 
see Eqs.|26l|3l|42]and|39] 



"Ac is the critical wave length of disturbances, 
''concerns the magnetic potential. 
''Aa is the grid spacing, 3D-disc. 

''axisymmetric limit, finite size disc (inner edge a,„, outer edge flout), 
symmetry with respect to the mid-plane, finite size layer (thickness 
2h = H), explicit function of vertical stractification, local validity. 

Table 1. Various prescriptions for the softening length adopted 
in some simulations of self-gravitating gaseous discs (see also 
Fig.D. 



0,5 



:5 



0,4 



0,3 





II 










nj 




II 

CC 




~1/e -0.3679 










DC 

i 1 



-3 -2 -1 



Fig. 5. Softening length normalized to the local semi-thickness 
h versus x as computed from Eq.[39]for h/a = 0.1 in the homo- 
geneous case. 

for all R. Only a numerical approach can yield the exact value 
of A, if it exists. However, a good approximation to this root 
can be obtained by considering the expansion of the complete 
elliptic integral of the first kind over its complementary modulus, 
as considered in Sect[3] To the lowest order, Eq. l40lbecomes 



where x is given by Eq. 



32] In other words, this is 
mi — 4e^^ . 



(41) 



(42) 



We therefore conclude that the appropriate prescription for A, 
valid over a few disc thicknesses around the singularity R - a, 
is given by Eq. |39] where m'^ is found from Eq. |42] 
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0,368 



0,3675 



0,367 - 



: 0,3665 



0,366 



0,3655 L 
0,001 




E=h/2a 



Fig. 6. The softening length at the minimum x - Q for the verti- 
cally homogeneous profile. 

4.3. Results for various stratifications 

It follows that A is proportional to h, and depends both on x and e. 
It is also sensitive to the vertical stratification through the func- 
tion Figure |5] displays A/h versus x for h/a = 0.1 in the ho- 
mogenous case (i.e., for;^' = Xo^ or 6xq - 0). We see that, in the 
range of validity, A goes through a minimum for x = 0. There, 
we have A:^ x h/2a and /o ~ 1/e and then the minimum is 



1 



1 



1 



1 + 



V2£ 



(43) 



This value is almost insensitive to the disc aspect ratio, pro- 
vided that the disc is geometrically thin. This is shown in Figure 
|6] where we plot A/ hat the minimum versus e. We note however 
that the main variations of A occur over a radial range of the 
order of2h (i.e., the total disc thickness). The thinner the disc, 
the sharper the variation around the singularity. 

Figure Q displays A/h for the mass density profile given by 
Eq.l20land for q = {1, 2, 3, 5, oo). At x = 0, we have the general 
formula 

A ^ e-(i+*^o) 
h ~ 



(44) 



Since 6xo = ^^^^ case, we also have 



22+2 1 

= - X I 



(45) 



4.4. An example of vertical mass density profile expandable 
in infinite series 

To illustrate the generality and power of the result, we consider a 
cosine profile of the form p(z) = cos which is a typical exam- 
ple where matter distribution is expandable in infinite series of 
the altitude. Actually, this profile can also be expanded by means 
of Eq.[3T|with the following coefficients : 



(-1)" 
(2«)! 



(46) 



0.5 



%0A 
-0.3679 

0.3 
• 0.2636 

0.2 





1 




1 


00 
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\ / 5, 


(homogeneous) 

- 








^\\\\ / // q=1 (parabolic case) 










Cosine 
















^ R=a- 


\/ ^ 

¥ + 




-0.2539 

1 


R=a 





Fig. 7. Softening length normalized to the local semi-thickness 
h versus x as computed from Eq.|39]for h/a - 0.1. The density 
profile is according to Eq.|20]with q - {1, 2, 3, 5, oo). The cosine 
profile is also shown. 



Figs. |4] and Q for ;^f and Ajh respectively. In particular, at x = 0, 
we have 



{-lYq 



(Iq+miq) 



371... 



(47) 



which results in \ 

h 

be obtained from Eq. 
parabolic case). 



e"' X 0.690 ~ 0.2539 (this value would 
]for q a: 0.849, which is close to the 



5. Concluding remarks 

In this paper, we have reported the first reliable prescription for 
the softening length A to be used in the numerical determina- 
tion of the gravitational potential of geometrically thin discs 
within the framework of softened gravity. This expression has 
been found by rigorously comparing the Newtonian potential of 
a geometrically thin disc (of finite thickness) with the softened 
potential of a flat disc. This is a function of the radius and disc 
local aspect ratio, and also depends on vertical stratification. It is 
accurate at the singularity and around (typically a few disc thick- 
nesses in radius). Although this formula is valid only locally, and 
obtained in the axisymmetric limit, it should help to improve the 
quality (realism, accuracy, and computing time) of 2D- and 3D- 
simulationfl If necessary, the accuracy of the prescription can 
be improved by considering higher order terms in the expanded 
complete elliptic integrals of the first kind. 

Finally, it would then be interesting to generalize our re- 
sults i) to mass density profiles that extend to infinity, such as 
the Gaussian profile (which corresponds to vertically isothermal 
discs), and ii) to the entire disc since this prescription is expected 
to be inaccurate far from the singularity. 
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